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Q i Abstract 

We discuss the generalized eigenvalue problem for computing energies and matrix ele- 
ments in lattice gauge theory, including effective theories such as HQET. It is analyzed 
how the extracted effective energies and matrix elements converge when the time sepa- 
rations are made large. This suggests a particularly efficient application of the method 
for which we can prove that corrections vanish asymptotically as exp(— (_E7jv+i — E n ) t). 
The gap -E/v+i — E n can be made large by increasing the number N of interpolating 
fields in the correlation matrix. 

We also show how excited state matrix elements can be extracted such that contami- 
nations from all other states disappear exponentially in time. 

As a demonstration we present numerical results for the extraction of ground state and 
excited B-meson masses and decay constants in static approximation and to order l/m-b 
in HQET. 
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1 Introduction 



In the early days of lattice gauge theories, K. Wilson suggested to use a variational 
technique to compute energy levels in lattice gauge theory [1]. The general idea was 
quickly picked up and applied to the glueball spectrum [2, 3] and to the static quark 
potential(s) [4]. 

The usual choice of a variational basis starts from some interpolating fields Oi(xo) 
already projected to a definite momentum and other quantum numbers such as parity. 
(Examples are given in the section on numerical results.) The associated Hilbert space 
operators Oi define states 

|<fc>=6i|0> and |&> = e -*> A / 2 |&), (1.1) 

where to is some time parameter and H the Hamilton operator. A variational principle 
is formulated by (t > to) 

Clearly Xi(t, to) yields a variational estimate of the lowest energy, Ei, with quantum 
numbers of Oi via Xi(t, to) ~ e - ^* - * ). At the same time Xi(t,to) is the largest 
eigenvalue of a generalized eigenvalue problem (GEVP), which we define below. 

The GEVP is applicable beyond the computation of the ground state energy and 
it has been widely used in lattice field theories [5]. Recent examples are found in [6-18]. 
Corrections to the true energy levels decrease exponentially for large time t, when the 
energies are extracted from the generalized eigenvalues as described below [19]. 

Apart from the paper of Liischer and Wolff [19], statements about corrections due 
to other energy levels seem to be absent in the literature. We here add such statements 
and suggest a special use of the GEVP, which we will show to be more efficient under 
certain conditions. We also treat the case of an effective theory and show numerical 
results for HQET. Furthermore, we show that the eigenvectors of the GEVP allow to 
construct operators (see eq. (2.16)) which applied to the vacuum yield states converging 
exponentially to eigenstates of the Hamiltonian. 

With a proper condition for t, to we prove that the convergence rate is ~ e -( E N+i-E n ) 
for the energies and ~ e~( E i*+i— E n)to f or m atrix elements of the constructed states. The 
occurrence of the large gap E^+i — E n is the main point of this paper. An earlier 
presentation of our work on this topic can be found in [20]. 

An indication of the effectiveness of the GEVP in suppressing contaminating con- 
tributions from excited states can be seen in figure 1, where we plot the effective energy 
of the ground state of a static- light pseudoscalar meson as obtained from different levels 
of Gaussian smearing and from the GEVP; for details, we refer to section 4. 
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Figure l: Plot of aEf (t, to = 4a) against t/a as obtained from the GEVP for a 3x3 and 
a 4 x 4 system, and from the effective mass plot of two correlators with different amounts 
of Gaussian smearing (22 and 90 iterations, respectively). The improved convergence 
of the GEVP solutions can be seen clearly. 



2 The Generalized Eigenvalue Problem 
2.1 The basic idea 

We start from a matrix of Euclidean space correlation functions 

oo 

Cij(t) = (O i (t)O*(0)) = Y,e- Ent i>mr nj , i,j = l,...,N 

n=l 

4>ni = C0n)i = {0\Oi\n) E n < E n+ i . 



(2.1) 



Note that we have assumed non-degenerate energy levels. Space-time could be contin- 
uous, but in practice applications will be for discretized field theories. In that case, 
we assume that either the theory has a hermitian, positive transfer matrix as for the 
standard Wilson gauge theory [21,22], or we consider t,to large enough and the lattice 
spacing, a, small enough such that the correlation functions are well represented by a 
spectral representation and complex energy contributions are irrelevant [23,24]. States 
|n) with {m\n) = 5 mn are eigenstates of the Hamiltonian (logarithm of the transfer 
matrix) and all energies, E n , have the vacuum energy subtracted. Oj(t) are fields on 
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a timeslice t that correspond to Hilbert space operators Oj whose quantum numbers 
such as parity and flavor numbers are then also carried by the states |n) , n > 0. In the 
above formula, contributions due to a finite time extent of space-time are neglected, so 
the time extent T has to be large. 1 

Before discussing the GEVP, let us comment briefly on a wider scope of applications. 
Besides the energy levels E n one may want to determine matrix elements (0\P\n) or 
(m\P\n) of operators P that may or may not be in the set of operators {Oi}. We will 
see that the GEVP provides a systematic approach to this problem as well. 

The GEVP is defined by, 

C(t)v n (t,t ) = \ n {t,t )C(t )v n {t,t ), n = l,...,N, t>t . (2.2) 

Liischer and Wolff showed that one can use it to determine systematically also excited 
states, namely [19] 

E n = lim Ef{t,to), (2.3) 
Ef(t,to) = -d t logA„(Mo) = ~[log\ n (t + a,t ) -logA„(Mo)]- (2.4) 

For a while (until Sect. 2.2) we now assume a simplified case, which helps to under- 
stand the usefulness of E n s and v n . The simplification is that only N states contribute, 

N 

C ij {t)=c\f{t) = Y J Z- Ent ^nir nj - (2-5) 

n=l 

We introduce the dual (time- independent) vectors u n , defined by 

(u n ,ip m ) =S mn , m,n<N (2.6) 
with (u n ,tp m ) = YA=i( u *n)i^mi- Inserting into eq. (2.5) gives 

C<°>(tK = e- £ «Vn , C(°)(t) u n = \£\t, t ) C(°)(t ) u n . (2.7) 
So the GEVP is solved by 

Af(Mo) = e-^*-*") , v n (t,t ) oc u n (2.8) 
and there is an orthogonality for all t of the form 

(n m , C(°) (t) un) = 5 mn p n {t) , p n (t) = . (2.9) 

These equations mean that the hermitean conjugates of the operators 

N 

An = J2«)iOi = (u n ,d), (2.10) 
i=l 

1 See [25] for a discussion of finite T when the O t carry vacuum quantum numbers. 
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create the eigenstates of the Hamilton operator, 

|n) = it |0) , H\n) = E n \n) . (2.11) 
Consequently matrix elements are easily obtained. For example 

(0|P|n> = (0|Pit |0) , (2.12) 
translates into correlation functions of P, 0{ . 

2.2 General statements about energies and matrix elements 

Let us now come back to the general case eq. (2.1). The idea is to solve the GEVP, 
eq. (2.2), "at large time" where the contribution of states n > N is small. We then 
consider the effective energies, eq. (2.4), and effective operators 2 

if (Mo) = e-^Qf(Mo), (2.15) 
Qf(Mo) = Rn(6,v n (t,t )), (2.16) 

R n = MMo), C(tK(Mo)r 1/2 ¥^TTY • (2-17) 
Corrections to the large time asymptotics are parameterized by 

Ef(t,t ) = E n +e n {t,t ), (2.18) 

oo 

e-^(Qf(Mo)) f |0) = |n) + X>nn'(Mo)|n'}. (2.19) 

n'=l 

The terms e n , ix nn i will disappear exponentially at large times. Note that in the literature 
the energy levels are sometimes not extracted like that and rather the standard effective 
masses of correlators made from A n = (O,v n (ti,to)) are used with fixed t\,to and the 
question of the size of the corrections is left open. However, the form eq. (2.18) has a 
theoretical advantage as it was shown in [19] that (at fixed to) 

e„(t,to) = 0(e- A£ "*) , AE n = min \E m - E n \ . (2.20) 

This is non-trivial as it allows to obtain the excited levels with corrections which vanish 
in the limit of large t, keeping to fixed. However, it appears from this formula that the 
corrections can be large when there is an energy level close to the desired one. This is 



2 For simplicity we discuss a special case of a more general class of operators depending on several 
times, which then all have to be taken large. An example is 

Sf = R n (6,v„(ti,to)), (2.13) 

R n = (Mti.to), C(fe)"n(«i,«o))- 1/a A w°* fa / 2 :y ■ (2-14) 

A n (to + 12, to) 
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the case in interesting phenomena such as string breaking [26-30], where in numerical 
applications the corrections appeared to be very small despite the formula above. Also 
in static-light systems the gaps are typically only around AE n « 400 MeV, and in full 
QCD with light quarks, a small gap AE n ps Im^ appears. 

Our contribution to the issue is a more complete discussion of the correction e n to 
E n as well as the definition of the effective operator, eq. (2.16), and a discussion of the 
corrections ir nn r to its matrix elements. It turns out that in all this, a very useful case 
is to consider the situation 

t >t/2, (2.21) 

e.g. with t — to = const, or 2 > t/to = const, and take to ( m practice moderately) large. 
In this region of t, to the second order perturbations in higher states n > N are not 
larger than the first order ones. We can then show that 

e n (Mo) = 0(e- A ^+ 1 -"*), AE m>n = E m -E n , (2.22) 
tw(Mo) = 0(e~ AElf + 1 - nt °) , at fixed t- 1 (2.23) 
7ri n (Mo) = 0(e- AS ™*°e- AB2 < l( '-' o) ) + 0(e- AEN + 1 ' lt ). (2.24) 

The large gaps A_Ejv+i,n can solve the problem of close-by levels for example in the 
string-breaking situation, but also speed up the general convergence very much. For 
example in static- light systems AE^i ~ 2 GeV translates into a gain of a factor of 5 in 
time separation. 

Eq. (2.23) means that matrix elements of (time-)local operators P and the associ- 
ated field P can be computed via 

{0\&*e- At Pe- At (Qff?\0) = (Qf (2t)P(t)(Qff (0))*) = (n\P\n') + (e™ A ^+ 1 ."'°) 

(2.25) 

and similarly for (n\P\0). 

We now turn to an outline of the proof of these statements, delegating the technical 
part to an appendix. 

3 Perturbation theory 

We start from the solutions above for C = and treat the higher states as perturba- 
tions. This perturbative evaluation was already set up by F. Niedermayer and P. Weisz 
a while ago [31] but never published. We noted the role of to > t/2, the form of the 
corrections to the effective operators defined above and could show that these relations 
hold to all orders in the perturbative expansion. 

We want to obtain A n and v n in a perturbation theory in e, where 

Av n = X n Bv n , A = A<®+eAW, B = B<® + efl^ . (3.1) 
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All these N x N matrices are assumed to be hermitian, B and B^ are assumed to be 
positive. We will set 



= C(°)(t), eA« =C«(i), 
= C(°)(to), e J B« = C( 1 )(t ) 



(3.2) 
(3.3) 



in the end with 

oo 

C S } W = E e_£n Vn^-, (3-4) 
n=AT+l 

which is hermitian and positive under the assumptions listed in the previous section. A 
possible choice for e is e = e -( E N+i-E N )t which becomes arbitrarily small for large to- 
The solutions of the lowest-order equation 

A<%W = WbWv®, (3.5) 

satisfy an orthogonality relation 

(vM,B(%M)=p n 8 nm (3.6) 

as in eq. (2.9) above. Non-degenerate eigenvalues A^ > A^L are assumed. Expanding 

A n = A(°) + £ A« + e 2 A( 2 ) . . . , v n = vW + «,« + . . . (3.7) 

with the orthogonality condition 

( V W,B(%( >) = 0, fc = l,2,..., (3.8) 

we get for the first two orders 

A<®vW+AW V W = 4°) [BCJ^ + BW^)] +A« j B(%(°), (3.9) 

+A( l 2 ) J B(°) V (°). (3.10) 

With the orthogonality, eq. (3.6), of the lowest order vectors, one obtains just like 
in ordinary quantum mechanics perturbation theory the solutions for eigenvalues and 
the eigenvectors 



A« 
Ai 2 ) 



= Vc^ 1 

/ ., "to 



) w (°) a M = p- 1 

nm > "nm rm 



rn=£n 



v {0) A v (0)> ! 

j ^n ^n I 

7(o) 7(o) 

An Am 



(3.11) 
(3.12) 



(^,A^ 0) ) 



x(0) _ \ (0) 



Pn 2 (4 0) ,A„^)) («<?>, B«i;<°>) .(3.13) 
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A recursion formula for the higher-order coefficients is given in App. A. 

We note that for our case of interest, to can be chosen large enough such that the 
perturbation theory is absolutely convergent. It therefore suffices to discuss the large 
to,t asymptotics order by order in the expansion in order to show the detailed form of 
corrections, in particular that they are controlled by the large gap AE^ + i tn . Conver- 
gence of the expansion is guaranteed for the following reason. As in quantum mechanics, 
the perturbation theory can be set up with the help of the resolvent depending on a 
complex variable z. The positions of its poles yield the (here generalized) eigenvalues. 
The resolvent for the GEVP is 

oo 

G(z) = (zB - A)" 1 = G (z) £ c k ~ zBW) G (z)] k , (3.14) 

k=0 

in terms of the unperturbed resolvent Gq(z) = (zB^ — A^)" 1 . The series expansion 
eq. (3.14) converges for large to (and t > to) since \\e(zB^ — A^)\\ < \\Gq 1 (z)\\ holds 
at sufficiently large to (note that the resolvent is needed only away from the zeros of 
G \z)). 

3.1 Application to the perturbations C^> 

Now we insert our specific problem eq. (3.2), eq. (3.3). With a representation (for 
m > n) 

( A (0)_ A (0) r l = {x (0) r l {1 _ e -(E m -E n )(t-t ) r l 

oo 

= (A^)- 1 ^-^-^*-* ), (3.15) 

1=0 

and some algebra spelled out in App. A. 2 one can show that eqs.(2.22 - 2.24) hold at 
any order in e. 

In particular, we list below the explicit expressions up to second order in e, keeping 
only leading terms in the sums. They demonstrate how the condition to > t/2 comes 
about. 

For the effective energies we have 
e n (Mo) ~ e-V»+ 1 -^ t [l-e-( E »+ 1 - E ») a ]c n , n>N+1 (3.16) 

_|_ ^ ^ g— (2Ejv+i— E n — .E m )tOg— {E m — E n ){t— to) _ g— (E m — E n )a 
m>n 

E —(2E N+ i—E n —E m )to —(E n —E m )(t—to) -i _ —(E n —E m )a i „ _ I 2 

c c 1 e | n,m,jv+l I 

m<n 



I |2 
\C n ,m,N+l\ 



where the first term (line) comes from order e in the perturbative expansion and the 
second and third terms come from order e 2 . Also, we have defined 

c n ,m,l = (u n ,ipi){ipi,u m ) . (3.17) 
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We see that, due to cancellations of t- independent terms in the effective energy, the first- 
order correction is independent of to and exponentially suppressed in t with the large 
energy gap En+i — E n as coefficient. (This correction is positive, while the sign of the 
combined second-order one depends on the energy shifts for the states around E n .) At 
fixed to and asymptotically large t the second-order corrections dominate because they 
follow the slower decay eq. (2.20). However, simple inspection of the above corrections 
shows that the to-dependent prefactors supress the second order correction sufficiently 
when to > t/2. Then the first-order one dominates, confirming eq. (2.22). 

For the amplitudes 7r nn '(t,to) we have main contributions at first order in e for 
n' n , n' < N. These contributions from n' > n and n' < n are given respectively by 

vw(Mo) = -V,iv + ie- (£N+1 - £ " )io e- (£ »'- £ " Kt - to) (3.18) 

and 

tw(Mo) = c n , n ,, N+1 e-( EN ^- E "^ . (3.19) 

Note that in the n' > n case the value n' = n + 1 has a leading exponential behavior, 
whereas for n' < n all values of n' have the same exponential behavior. 



3.2 Effective theory to first order 

In an effective theory such as HQET, all correlation functions 

Cy(t) = qft(i) +coCl /m (t) + 0(u, 2 ) (3.20) 

are computed in an expansion in a small parameter, u, which we here consider to first 
order only. The notation is taken from HQET where to oc It also helps to avoid 

a confusion with the perturbation theory in terms of treated before. 



3.2.1 Energies 

We start our analysis from the GEVP in the full theory, eq. (2.2), and use the form of 
the correction terms of the effective energies it < 2io) 

Ef (Mo) = a" 1 log . y_f'* o) , . = E n + 0(e- AE ^), (3.21) 
A n {t + a, to) 

see the discussion above. Expanding this equation in u>, we arrive at 

E?(Mo) = K"-""{t,t< 1 ) + ^K" A, " , <t,t«) + 0(u,' 2 ) (3.22) 

Ef^M = »-Mo g j|ML (3 . 23) 
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with the behavior at large time, 

Ef ' stat (t, t ) = E% at + PT' e- AE %+i>" * + ..., (3.25) 
Ef^(t,t ) = + [/?V m -rtA£i / ; ii Je- A ^» t + .... (3.26) 

Following the beginning of Sect. 3 we now take the lowest order correlator to define the 
unperturbed GEVP, 

c stat (tK tat (Mo) = A^ at (t,to)C stat (to)< tat (^io), (3.27) 

whose eigenvectors with normalization (w^ at (t, to) , C stat (to) v^^t, to)) = 5 mn are then 
needed in the formula 

X^t {t ] t ° l = (< tat (^*o), [[Ar(i,io)]- 1 C 1 / m (t)-C 1 / m (io)]< tat (t,*o))(3.28) 
for the first-order corrections in u. We note that the finite time corrections in eq. (3.24) 

^^stat £ 

comprise of both a term e jv+i,™ anc j one w jth an extra factor of t, however only 
their prefactors are new, once AEfj-^ n is known from an analysis of the lowest order 
eigenvalues. Due to the large gap AEfj-^ n , the corrections disappear quickly with time. 

3.2.2 Matrix elements 

The discussion of matrix elements starts from the operators Q^ ff , eq. (2.16). Their 
expansion 

Qf = RT (vT\t, t ) , 6) +u [/# m « tat (t, t ) , 6) + RT («y m (t, to) , 6)] (3.29) 
is given in terms of 

RT = « tat (Mo), C(t),r(t,t ))' 1/2 X l^l'^ , (3-30) 

Rl/ m l K tat (*,*0),g 1 / m (tX at (t,t )) 

i^ tat 2 (^ tat (t,to), C stat (*X tat (Mo)) 

An //m (tp + t/2, t ) _ Ay m (t + t, t ) 
A^ at (*o + t/2,t ) A^ tat (io + Mo) 

and the first order perturbation to v n , 



N 

" A-t(Mo)-Ar(Mo) • (3 ' 32) 
A simple matrix element is then 

pf(t,t o ) = {0\Qfe-" t P\0) = R n (v n (t,t ),C P ), (C P )j = (O,-(0)P(t)> .(3.33) 



(tf at , [C 1 /-^) - A^Mo)*? 1 /™^)] < 
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For example when P denotes the time component of an axial current, it yields the decay 
constant of the associated ground state meson via, 



with the w-expansion 
pf(Mo) = 
K ff ' 1/m (Mo) = 



eft ,stat i 
Pn 



Pi = (0|P|1), 



(t,t )(l+u;pf' 1 / m (t,to) + 0(^)) 



R. 



l/m 



« at (Mo),4 /m (t)) K /m (Mo),Cff at (t)) 
R*°* ' (v^(t,t ),C^\t)) + (v^(t,t ),C^(t)) 



l/m , 



+ 



(3.34) 



(3.35) 
. (3.36) 



Specializing further, the large time asymptotics of the effective ground state matrix 
element has the following form 



pf(Mo) ~ Pl+1 e- AE »+^ t0 



Pi (Mo) 



p stat + y 



eff, l/m 
Pi 



(Mo) 



stat e -AEf + ' u i 

p V m + e-^fi.i*" 



y/m _ 7 stat p l/m 



to A^ +lil 7 



(3.37) 
(3.38) 

(3.39) 



We note that the correction term proportional to toe n+i,i ° [ s entirely fixed by a 
first order analysis of the energies and the lowest order one of the matrix element. 



4 Application to static-light B s -mesons 

We have carried out numerical tests of eqns. (3.23)-(3.39) in the context of static-light B s 
mesons in quenched HQET with a HYP2 static quark [32-34] and a non-perturbatively 
0(a)-improved [35,36] Wilson valence (strange) quark 3 . Our lattices are L 3 xT with 
L ~ 1.5 fm, T = 2L and periodic boundary conditions in all directions. For this 
demonstration we have chosen a = 0.07 fm {(3 = 6.2885, k s = 0.1349798) and use an 
ensemble of 100 quenched configurations. A more detailed analysis including multiple 
lattice spacings to enable us to take the continuum limit is in progress. 

Strange quark propagators are computed using a variant of the Dublin method [37]. 
We use approximate instead of exact low modes and employ even-odd preconditioning 
in order to reduce both the size of the eigenvalue problem to be solved for the low modes 
and the number of noise sources to be used for the stochastic estimator by a factor of 2; 

3 One may object that a quenched calculation does not satisfy our basic premise, eq. (2.1). However 
the quark considered here is rather heavy and experience has shown that quenching represents a small 
modification of the full theory in such cases. Clearly a demonstration in a unitary theory will be 
welcome. 
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this leads to a significant reduction in effort and is expected to also reduce the amount 
of noise introduced by the stochastic estimation of the short-distance contributions [38] . 
The interpolating fields are constructed using quark bilinears 

O k {x) = Mxhol^l k \x) (4.1) 
Otix) = ? 1 (fc) (x) 70 75V'h(^) (4.2) 

of identical Dirac structure, but with different levels of Gaussian smearing [39] for the 
light quark fields 

4 k) (x) = (1 + k g A) r * U*) , (4-3) 

where the gauge fields in the covariant Laplacian are first smeared with 3 iterations 
of (spatial) APE smearing [40,41]. We use R k = 22, 45, 67, 90, 135, 180, 225 with 
kg = 0.1. The local (Rq = 0) operator is also included in order to be able to compute 
the decay constant. Here and throughout, iph(x) denotes the static quark field. 
For these operators, we compute the following correlators: 

„13 

cf\t) = —J2(oMt + t x )o*(x,t x )), 

Cf A {t) = (Oi(y,t + t x )OUx,t x )), 

S,y,t x 

C]f(i) = £ (O t (z,t + t x )O kin (y,t y )0*(x,t x )) , 

C*f n (t) = £ (Oi(z,t + t x )O spin (y,t y )0*(x,t x )) , (4.4) 

where 

SA {x,t) = Vh(^>*)7o75[7 • Dip[ 0) ](x,t), 

O kin {x,t) = i/j h (x,t)[D 2 ip h ](x,t), 

O spin (x,t) = ^ h (x,t)[a ■ Bip h ](x,t). (4.5) 

The correlator C SA is introduced to ensure the axial local current is correct up to both 
O(a) and 0(l/mt,), whereas the last two correlators are the C(l/mb) terms from the 
HQET action. 

The energy of a state \n) is given by 

E n = ET l + Wki„^ n + io spin E s T r + Sm (4.6) 
and its decay constant factorizes as 
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where Z^ ET , ct>kiiu Wspin> c^ ET and 5m are (divergent) matching constants between 
the QCD action and currents and their HQET counterpart. A strategy to compute them 
has been presented in [42,43] and we quote their approximate numerical values [38] at 
= 6.2885 (Z" QET ~ 1, tJkin/a ~ 0.4, w spin /a ~ 0.7, c" QET /a ~ -0.6) only as an 
illustration of the expected size of their contributions. 

The resulting 8x8 correlator matrices C(t) are symmetrized and then truncated 
to iV x N matrices C ( - NxN \t) by projecting with the N eigenvectors belonging to the 
N largest eigenvalues of C stat (tj): 



C Stat (ti)bn 



X n b n 

biC{t)b m , n,m<N. 



(4.8) 
(4.9) 



For N not too large, this helps to avoid numerical instabilities in the GEVP that could 
otherwise lead to large errors [44]. We use tj = 2a. 

For each of the resulting N x N correlators, we solve the static GEVP and com- 
pute the static and 0(l/mb) energies and matrix elements as per eqns. (3.23)-(3.24) 
and (3.30)-(3.36). This gives a series of estimates £^' stat (t, to), p^' stat (t, to) etc. with 
associated statistical errors, which we determine by a full Jackknife analysis. 

To arrive at final numbers for E n , p n , we first need to estimate the size of the 
systematic errors coming from the higher excited states. To do this, we perform a fit of 
the form 

pstat ^stat ^£ 



(4.10) 



(cf. eq. (3.25)) to the GEVP results for En ,sUt (t,t ), fitting the data at TV = 3, . . . , 5, 
to/a = 3, ... ,6 and n = 1, . . . , 6 simultaneously. Then, using the values of £^ tat and /3^ ta ^ 
determined from this fit as (fixed) input parameters, we fit En' km {t, to) and £^' spm (t, to) 
by 



E^' kin (t,to) 



e -{E$* 1 -Eg>*)t ( 4 _ n > 



£^ spin (Mo) = Et n + 



(E^-E^t 



/okin astat + / 771km 771km 
Pn,N ~ Pn,N 1 {-^N+l ~ 

(cf. eq. (3.26)) in the same manner. Subsequently, we also fit p^' stat (t, to) to 

p ^\t, t ) = P T + 7 s>-( E " a +i- £ « at )'° 

(cf. eq. (3.38)), and j>n' km (t, t ), pn ' spm (t, to) and Pn' 5A (t, to) 



(4.12) 



(4.13) 



to 



Pn'^Mo) 



JV,spin/ 



„,kin stat kin „,stat j. / 771km zrikim 

ln,N ~ ln,NPn ~ ln,N *0 [E N+1 ~ E n ) 



(4.14) 



p^ SA (t,to) = a a + 



SA 



SA stat <5A 
ln,N ~ ln,NPn 



-(B|tat i _ £ Btat )t 



(4.15) 
(4.16) 
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t 

Figure 2: Plot of £^' stat (t, to = 4a) against t. All quantities are given in lattice units. 
Data points are displaced horizontally for better visibility. Also shown for the purposes 
of comparison are the fit curves from eq. (3.25) and the plateau given in table 1. Note 
that the curves are fitted to the total data set, not just the points shown in the plot, 
and that the plateau is obtained at a different to than shown here. 

(cf. eq. (3.39)), using the previously determined values for E^ at , 7^y, E^ in and E s ^ m 
as fixed input parameters. 

Since our fits do not include contributions beyond the leading excited state cor- 
rection, the errors estimated from the fits may not be very reliable. In particular we 
consider the fitted values of £'| tat , . . . , £ , | tat as rough estimates only, which may have 
significant additional systematic errors. We therefore use the fitted values only for the 
purpose of determining the systematic correction e(t) coming from the excited states. 
The good quality of the fit that can be seen from figs. 2 to 9 makes them very suitable 
for that purpose. In particular we note that we find the E^ (t, to) to be essentially 
independent of to, as predicted by our fit formulas. 

For a reliable estimate of our quantities of interest, we calculate plateau averages of 
the energies from t = t m i n > to to t = 2to at each N and to (and from to = *o,min >t/2 
to to = t at each N and t for the matrix elements). As our final estimate we take that 
plateau average for which the absolute sum <7 tot = Cstat + ^sys of the statistical error 
o"stat of the plateau average and the maximum systematic error fj sys = e(t m j n ) becomes 
minimal, subject to the constraint that cr sys < cr s tat- We impose the latter constraint in 
order to ensure that the systematic errors are subdominant. 

Figs. 2 to 9 show plots of our results at different values of N along with the fits 
and the optimal plateau and its error bands. It is apparent not only that the data are 
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Figure 3: Plot of E 2 ,stat (t, io = 4a) against t in the same style as fig. 
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Figure 4: Plot of E 1 ' ]aD (t, to = 4a) against t in the same style as fig 
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Figure 5: Plot of E^' spm (t,to = 4a) against t in the same style as fig 
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Figure 6: Plot of p l ' stat (io + 3a, to) against to in the same style as fig 
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Figure 8: Plot of p^' km (to + 3a, to) against to i n the same style as fig 
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Quantity 


run nx 


Plateau fit 


Plateau 


Precision 


„ ZT'stat 


0.3U2o(4J 


U.o0o6(yj(z J 


5,9,12 


P^o Mev 




0.483(3) 


0.484(7) (5) 


5,6,10 


~ 20 MeV 


aEf at 


0.61(1) 


0.59(3) (2) 


5,6,11 


~ 100 MeV 


aEf at 


0.66(1) 








a£f tat 


0.74(2) 








aE stzt 


0.79(4) 








a 3/2 p stat 


0.0689(2) 


0.0694(5)(2) 


5,8,12 


0.7% 


a 3/2pStat 


0.086(3) 


0.089(10)(4) 


5,7,10 


12% 


a 2 Ef n 


0.6731(4) 


0.6729(4)(1) 


3,5,6 


<0.6 MeV 


a 2 Ef n 


0.691(5) 


0.69(2)(2) 


4,8,11 


~ 40 MeV 


a 2 Ef n 


-0.0129(1) 


-0.0131(3)(1) 


4,8,9 


<0.7 MeV 


a 2 Ef n 


-0.0106(4) 


-0.0113(9)(6) 


5,6,7 


~ 2 MeV 


apf n 


-0.622(8) 


-0.619(12)(3) 


5,5,8 


0.6% 


spin 

aPi 


0.385(2) 


0.382(3)(1) 


5,6,8 


0.7% 


ap{ A 


0.3120(3) 


0.315(3)(1) 


4,6,18 


0.1% 



Table l: The values of quantities of interest as determined from the fit and plateau 
picking procedure described in the text. The errors in the column headed "Plateau fit" 
are statistical and systematic. The column headed "Plateau" contains N, to and t nun 
of the best plateau for energy levels, and N, t niin and t of the best plateau for matrix 
elements. The column headed "Precision" contains an estimate of the contribution of 
the total error of that quantity to the total absolute error on the energy of the state, or 
to the total relative error of the decay constant, respectively. Note that terms such as 
<^kinPi m do not directly give the physical 1/m correction since they contribute divergent 
contributions cancelled for example by Z^ ET . 
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Figure 9: Plot of Pi ,spm (to + 3a, to) against to in the same style as fig. 2. 

very well described by the fitted functional forms, but also that the fits agree within 
errors with the more conservative plateau estimates. 

To study how independent the results obtained from the GEVP are from the op- 
erator basis used, we have rerun our analysis for the static case using a different basis 
of operators where we exchanged some of the intermediate smearing levels against op- 
erators constructed by an additional application of the Laplacian. We found that the 
results for both the static energies and the static matrix elements were the same within 
errors as those obtained from our original basis. 

5 Discussion 

We have performed a theoretical analysis of the GEVP for lattice field theories, eq. (2.2). 
The N x N problem is expanded in terms of a convergent perturbative expansion around 
the unperturbed system defined by truncating the spectral representation of the cor- 
relation matrices after N levels. The contribution of all the higher levels defines the 
perturbation. 

The GEVP involves two time separations, to,t, where t > to- At large t and fixed 
to, the extracted energies E^ 1 converge to the eigenvalues of the Hamiltonian with a 
rate 

0(e~ AEnt ), (5.1) 

where AE n is the distance of E n to the closest energy level. Such a spectral gap is 
typically only a few hundred MeV, requiring t to be above t > 1 fm for an appreciable 
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suppression. One realizes, however, that the first order corrections are only of size 




(5.2) 



where En+i — E n can be much bigger than AE n . At the second order in perturba- 
tion theory and beyond, the mixing of different levels introduces the slower decaying 
corrections, eq. (5.1) when to is kept fixed, t is taken large. We could show that the 
favorable suppression, eq. (5.2), is recovered to all orders if one chooses to > t/2. This 
is possible because the mentioned mixing is suppressed exponentially in to- Thus, even 
if it is challenging numerically, a large to (and N) is more important than a large t to 
keep systematic errors small. In the numerical demonstration, section 4, we have shown 
that this property can be used in practice. 

Furthermore, because of the exponential suppression of mixing with to, one can also 
write down interpolating fields, eq. (2.16), for all states with energies up to E n = E^. 
These fields are linear combinations of the N fields one started from, with coefficients 
which depend on to- They are interpolating fields for the desired energy eigenstates up 
to corrections which decay as 



again with the large energy gap. 

In a numerical example we could demonstrate eq. (5.2) and eq. (5.3) nicely by 
increasing N, to,t starting from small values. The decrease of the correction terms is 
clearly observed, in particular in Fig. 2 and Fig. 6, and the fit demonstrates that the 
derived formulae are (at least approximately) valid at the accessible values of N,to,t. 
This quantitative understanding of correction terms allows for a determination of energy 
levels and matrix elements with subdominant systematic errors. 

In the considered HQET problem we achieved a sub-percent level (statistical + sys- 
tematic) determination of the matrix elements needed for the B- meson decay constant, 
both at the leading and at the next-to- leading order in 1/mb- The ground state energy 
is obtained at the level of about 3MeV. A more complete analysis at several values of 
(3, including a continuum extrapolation, is currently in progress. 
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A Recursive perturbative expansion 
A.l The general case 

Solving eq. (3.1) order by order in e and considering the result at order e k , the GEVP 
reads 

= (A<°>-A£>B<°>)t,<*> + (A^-Al 1 )^))^- 1 ' 

+ (-X^bM - AP)B(°>)t4 fc - 2 ) + (-Ai 2 )BW - Ai^C")^- 3 ) (A.l) 
+ ... + (-X^B^-X^B^, 

where 

A^^-Ai ^ 1 '. (A.2) 

Projecting with gives 

Al fc) Pn = (v^ , A n v^) - £ A« (t,<°> , flW^- 1 -") , (A.3) 



and a projection with , m ^ n yields 



- £a<M> (A.4) 



m^n 



a (fc) fl _ ^(0) B (0) (fc)N 



fc-1 



- X m { i =1 



z=i ; 

The combined recursions 

a (*) - o" 1 !( v (0) A -.(k-ih _ VA (Z) ^t; (0) B^ti^-'^ + fl a (fc_/) 

"mn — Hm ( ) (o) ) V m ' ^"^n J Z—t n ^ m ' n ) ' Hm^mn 

An — A m ; =1 J 

(A.6) 

and eq. (A.3) then determine the solution to arbitrary order in the perturbations starting 
from the initial values af!, = S mn , Xn\t,to) = e _E «( t_t °). 



A.2 The case of Euclidean QFT 

We now apply this with 

A® = c(°)(t), 

B<® = C(°)(t ), 
v {0) - u -> 



eA( 1 ) = C'W(t), (A.7) 
eB (i) = C«(t ), (A.8) 
/On = e"^* , (A.9) 
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where we recall that 



N 



= Ee-^W;,, C\]\t) = ^ (A.10) 

n=N+l 



n=l 



and 



^m = {n\di\0), H\n) = E n \n). 



(A.ll) 



This means that eAn ^ • • • are corrections to the desired quantities which are the 
lowest order u n ,)$\ We also recall the orthogonality 

(u n ,ip m ) = 5 mn ,m,n<N. (A.12) 
A.3 Proof of eqs.(2.22 2.24) to all orders 

The first step is to identify what exactly has to be shown concerning the behavior of 

(k) (k) 

An and q„ to prove eqs.(2.22 - 2.24). For e n we expand 



mmo) = -d t i og [i + ^ k ^wtr\ 

k>l *n (t, to) 



Ej^fl J ^n^jt, t ) \ [ n l) (t,to) 



fcl,...,fc;>l 



Ai 0) (Mo) AHMo) 



(0), 



.(A.13) 



If the conditions 



^A|J(Mo) = ^ e -A^ +1 , n t^ ^n ] (t,t ) 



An 3 '* (Mo) 



a! ^, to) 



0(1) 



(A.14) 



are satisfied, one sees easily that each term in eq. (A.13) is of order e AE N+i, n t_ j n 
other words eq. (A.14) is a sufficient condition for eq. (2.22) to hold. Next we expand 



af)t|0)=e-^(Qf)t|0) 



k>0 n'>l 



x [l + e n ] (A.15) 



1+ £ e fcl + fe2 Cfcl , fc2 (e)(t,to),e^*C(°)(t)^)(t,to)) 

fcl,fc2>l 



fel,fe2>l 
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with some irrelevant coefficients k 2 > Cfci ,fc 2 and with e n a shorthand for sums of terms 
£n(t', to), which are all negligible provided eq. (2.22) holds. We express the various terms 
in eq. (A. 15) through a mn in order to arrive at conditions for these coefficients: 



EE e ~ s "V>(vv,4 fc) (Mo)) 



k>0 n'>l 
N N 

= £ £ e " (£ "'" En) V> £(^,^)a^(t,to) + 0(e- A ^*) 



=1 m=l 



fc>0 



N 



= £e-(^-^|m) E«^(t,io) + O(e- A ^-' 



m=l 



fc>0 



A' 



= |n>+ £ e-^-^|m)E«^(^*o) + 0(e- AE ^*) 



(A.16) 



(^'.WWtflf 2 ') = 0(e- AEjv + 1 ."*) (A.17) 

(i;<fr\e^*C( > = if fci = or fc 2 = (A.18) 

m,m' 

= E«^ ) (*'*o)«S(t^o)e-^-^ t (A.19) 

TO 

A look at these terms shows that it is sufficient to prove 
aW(Mo) = 



m > n 



5 A^(Mo) 
*A^ 0) (i,t ) 
A^(Mo) 
Af(Mo) 



0(e 



0(1) 



(A.20) 
(A.21) 
(A.22) 



for all fc. (The last two conditions are sufficient to show that e n (Mo) = 0(e~ A - Ejv + 1 -™ *)). 
Note that for k = the above conditions hold trivially. 

As a next step we collect the large time behaviour of the different terms appearing 
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in the recursions eq. (A. 3) and eq. (A. 6) for A^ fc ) , a^: 

4°) = 0(1), (A.23) 
e£W = 0(e- EN+lto ), eA^ = 0(e- EN+lt ), (A.24) 
eA„ = o( e -EN+it e -E n (t-t )^ + ( e -^+i*), (A.25) 

e An _ Q^ e -(£ JV+1 -E m )t ^ _j_ Q^ e -(E JV+1 -B m )t e -(E N+1 -E n )(t-to)^ (A 27) 

where in A„ we keep the two terms because the leading one has no t-dependence. For 
shortness we drop the O symbol from now on but just count orders and we use t > to- 
This means e.g. that e~ EN+lt ° + e~ Ent ° = e~ Ent °. However, since derivatives with 
respect to t are relevant we have to be careful not to drop e~ Ent with respect to e~ Ent ° . 
The shorthands 

r)Nm(t) = 0(e-^+ 1 - B -)*) (A.28) 

g— E„t 

W*) = ^Ej—^Tt (A. 29) 



J 3 

are useful to discuss the large t, to asymptotics, 

AO) 



e - tJnt — e 

YljLo e~ j( - Em ~ En ^ when m > n 

~ YZLi e- j< - En ~ Em ^ when n > m 



(A.30) 



Xk> u [ l|0(e-( £ ™- £ ») f ) whenm>n 

lnm{t -to) = < nf~-(E n -E m )U ,„U™ „ ^ ™ ( A ' 31 ) 



Ai 0) - a£? ~ ~ I 0(e-( £ «-^)*) when n > m 



x (o) _ x (o) 

^n />m 
X (0) _ x(0) 



rjNm'{to) [1 + Wn(i - *o)] 7nm(* - *o) (A.32) 
mm' {to) lnm{t ~ to) ■ (A.33) 



In this notation the asymptotics of the first order corrections are 



(0) = r]Nn{to)[l+VNn{t-to)], (A.34) 

An 

£«mn = mm{to) [1 + VNn{t ~ *o)] 7nm(* ~ *o) , (A.35) 
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while the recursions are 

(k) >(i) fc-i x (0 

e (Q) - e Q) a +e ^ ( ) « 

x (i) r fc-i x (o ^ 

= ^{e k + e*" 1 £ ^ Wn (t ) a**" 1 " 1 ) (A.36) 

An [ 1=1 An J 

fc-1 ,(0 

e k a« = e fc a« a^" 1 ) + e*" 1 £ w m (t ) a**" 1 "') 7 nm(* - to) 

i=l An 

fc-i ,(/) 

+ efe E^)«^ 7n m (t-t0) 
Z=l An 

f fe_1 \W r i 1 

= a« a^- 1 ) + £ ^) W*o) + e k a^} (A.37) 

I J=l An J 

where a( fc ) = max m amn- We finally will need the start values of the derivatives 

ed t-Jy = VNn(to) T] N n(t - t ) = T] Nn (t) , (A.38) 
An 

e 9ta^ = rj Nm (t ) [da nm {t - t ) + d t r] N n{t - t ) j nrn (t - t ) 

+ VNn{t ~ t ) dt1 nm (t - t ) ] 

= mm^e-^-^-^, (A.39) 
where eq. (A. 31) has been used. In particular the relation 

eVNn(to)d t a$ n <r) Nn (t) if t >t/2 (A.40) 

or equivalently 

ed t a^ n <VNn(t-t ) if t >t/2, (A.41) 

will be relevant. 

In terms of the abreviations 

r k = e k ~ 1 ^, x k = e k a%l, X k = e k = e k max a W (A.42) 

An 

we want to show that 

Xfc = r] Nm (to) 7„m(t - *o) , r k = 0(1) , (A.43) 

d t x k = Wn(t-t ), d t r k = T] Nn (t - t ) , (A.44) 

which immediately implies eq. (A. 20) to eq. (A. 22) and thus proves eq. (2.22) to 
eq. (2.24). For k = 1, eq. (A.43) and eq. (A.44) are satisfied by the start values derived 



24 



above. The induction steps from k — 1 to k for k > 2 are very simple (remember we 
assume to > t/2): 



fc-i 



r fc = + r)Nn(to) = maxr?Ar m (t )7 nm (t - t ) = 0(1) , 

Z=l 

(A.45) 

a^fc = a;i |rfe + ^rix fe -i| = r)Nm(to) lnm{t - t ) , (A.46) 



fc-i 



d t r k = d t X k _ x +VNn{to)^2d t [riX k _i_ l ] = rj Nn (t - t ) , (A.47) 



l=i 

k-l 



d t x k = dt[xi r k ] + ^ d t [xi n x k _i] = r] Nn (t-t ). (A.48) 
i=i 

They conclude our proof. 
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